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Abstract. 
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The five lectures introduce to the description of entanglement in many-particle systems 
and review the ground-state entanglement features of standard solvable lattice models. 
This is done using a thermodynamic formulation in which the eigenvalue spectrum of a 
certain Hamiltonian determines the entanglement properties. The methods to obtain 
it are discussed and results, both analytical and numerical, for various cases including 
■ time evolution are presented. 
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Entanglement in many-particle quantum states has been a topic of intense research 
in recent years with applications in numerics and interesting links to statistical physics. 



Q\ ', It is therefore excellently suited for an advanced course in a summer school. The 
following notes correspond closely to five lectures given in July 2011 at the International 
Institute of Physics in Natal, Brazil. They are based on a recent review article [1], 
but the material has been properly adapted to the purpose. Thus they contain more 
introductory examples and certain topics are presented in more detail. On the other 
hand, new material from the last two years, as well as supplementary notes have 
been added. Throughout the notes, the style is lecture-like with itemized statements. 
References are only given in direct connection with the problem at hand, show a 
preference of own work and should not be regarded as an exhaustive list. Compared to 
the version handed out in Natal, additional figures have been included and some editing 
took place. 
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1. Background and basics 

In this section, I summarize the basic features of entangled states and reduced density 
matrices and illustrate them with examples. For further details, see e.g. the short review 
[2]- 

1.1. Introduction 

Entanglement is a notion which goes back to 1935 when it was introduced by Schrodinger 
in a series of three articles (in German, the German word is "Verschrankung" ) [3] . At the 
same time Einstein, Podolski and Rosen discussed their famous "Gedankenexperiment" , 
in which they considered two particles with fixed total momentum and relative distance. 
Nowadays this is usually formulated with two spins, and this is also where one encounters 
entanglement first. Entanglement has to do with the features of quantum states and 
the information contained in wave functions. For a long time, it was a topic discussed 
mostly in quantum optics and for systems with few degrees of freedom. 

In the last 25 years, however, it has seen a revival with input from very different 
areas, namely 

• the theory of black holes 

• the numerical investigation of quantum chains 

• the field of quantum information 

In these cases one always deals with large systems and many degrees of freedom. 
In entanglement investigations, one asks the following question: 

• given a total system in a certain quantum state 

• divide it (in space, or in Hilbert space) in two parts (bipartition) 

• how are the two parts coupled in ? 

This is more general than looking at, say, a two-point correlation function. And there 
is a general way to answer this question, namely one can bring into a standard 
form, which displays the coupling. This is the Schmidt decomposition which we will 
discuss in a moment. To obtain it in practice, one uses quantities which determine all 
properties of a subsystem, namely reduced density matrices (RDMs). They also contain 
the information on the entanglement and will be the basic tool throughout the lectures. 

The states we will study are the ground states of models which, on the one hand, are 
solvable and, on the other hand, have a physical significance, like tight-binding (hopping) 
models or spin chains. As in other contexts, they serve as points of orientation which 
allow to study the features of the problem and to develop a feeling and an overall picture. 
My own interest arose in connection with the DMRG, where the entanglement turned 
out to be crucial for the performance of the method. Entanglement continues to play a 
role also in other algorithms and their design, and in this respect it has quite practical 
implications. But in these lectures, we shall be concerned essentially with the theory. 
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1.2. Schmidt decomposition 

Consider a quantum system in state and divide it into two parts 1 and 2. Then one 
can write 

\^)^J2 A ^n\K)\K) (1) 

m,n 

where the |^^) and are orthonormal bases in the two Hilbert spaces. 
Note that one has a double sum and that the matrix A m ^ n is in general rectangular, 
since the dimensions of the Hilbert spaces can differ. Nevertheless one can obtain a 
diagonal form via the so-called singular-value decomposition 

A = UDV (2) 

where U is square and unitary, D diagonal and V rectangular with orthonormal rows. 
This gives 

|*>= E U m,nD n , n V^l)\^l) (3) 

m,n,k 

Combining |\^) with U and with V one obtains with A„ = D rhn 

i*> = E a «i^>i^> ( 4 ) 

n 

This is called the Schmidt decomposition (Schmidt f907) [4]. For the history see section 
1.7. It has the following features 

• Single sum, limited by the smaller Hilbert space 

• New orthonormal sets |$") in both parts 

• X] l^n| 2 — 1 if |^) is normalized 

• Entanglement encoded in the A„ 

• Limiting cases 

Ai = 1, A n = for n > 1: only one term, product state, no entanglement 
A n = A for all n: all terms equal weight, maximal entanglement 

This refers to a particular bipartition and one can investigate different partitions to 
obtain a complete picture. Some standard bipartitions for one-dimensional systems are 
shown in fig. 1. 

1 2 ^ - 2 1 2 




Figure 1. Bipartitions: Chain cut in two halves (left), ring cut in two halves (centre) 
and segment in an infinite chain (right). 
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1.3. Examples 

We give here examples for the Schmidt decomposition in three different systems. 

(a) Two spins one-half 

l^i) = product state (5) 

|¥ a ) = o|+)|+) + 6|+)|_) 

= |+) a\+) + b\— } product state (6) 

|* 3 ) = a|+)|+) + d\-)\-) entangled state (7) 
All these states are already in Schmidt form. However 
|*4> = a|+)|+) + 6|+)|->+c|-)|+> 

= 1+) [fl|+> + 6|->] (8) 
is entangled, but not in Schmidt form, because the two states in subsystem 2 are not 
orthogonal. 

(b) Two large spins [6] 

Consider the ferromagnetic spin one-half Heisenberg chain with N sites 

# = -J^s n s n+1 (9) 

n 

All eigenstates can be written as = \S, S z ) with total spin S and z-component S z . In 
the ground state all spins are parallel, S = N/2, and S z can be chosen. Choose S z = 
and divide the chain in two halves. Then one can use angular momentum addition as 
illustrated in fig. 2 




to obtain 

S/2 

\S,Q)= J2 c M \S/2,M) 1 \S/2,-M) 2 (10) 

M=-S/2 



with the Clebsch-Gordan coefficients 

S\ S\ 



^{2Sy. (S/2-M)\(S/2 + M)\ 
This is the Schmidt form for this state. Its features are 
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• Only (S+l) terms, while dimension of Hilbert space is 2 s 

• For large S, cm ~ exp(— 2M 2 /S), Gaussian 

• Analogous formulae for arbitrary S z 

• Special case S z = S, all spins in the z-direction. Then \S, S) is a product state 

\S,S) = \S/2,S/2) 1 \S/2,S/2) 2 (12) 

(c) Two coupled oscillators [7] 

Consider the Hamiltonian (m = h= 1) 

H = \{p\ + ulx\) + ^ipl + ulxl) + ^k( Xl - x 2 f (13) 

The eigenfrequencies are uf = Uq + 2k and io\ = Uq with corresponding normal 
coordinates 

1 1 

Vl = ~ y2 = ^(x 1 + x 2 ) (14) 

In these coordinates, the ground state is the product of two Gaussians 

l*o> = (^) 1/4 ex 
Then the following formula holds 

-tanh 
cosh 77 

n=0 1 

where exp(4^) = uji/u 2 and the |<&„) are oscillator states for a frequency u = ^Uiu 2 , 
i.e. in between u\ and u 2 . 

Features 

• Schmidt states are "squeezed" states 

• Coefficients decay exponentially, \\ ~ exp(— en) 

• weak coupling k: U\ ~ oo 2 — > rj small, e large, rapid decay, weak entanglement 

• strong coupling k: oj\ ^> u 2 — > r\ large, e small, slow decay, strong entanglement 

These features are also found for one oscillator in a whole assembly. 
1.4. Reduced density matrices 

The Schmidt structure just discussed can be found from the density matrices associated 
with the state This is also the standard way to obtain it. Starting from the total 
density matrix 

p=\*)m (it) 

one can, for a chosen division, take the trace over the degrees of freedom in one part of 
the system. This gives the reduced density matrix for the other part, i.e. 

Pi = tr 2 (p) , p 2 = tri(p) (18) 



*o> = ^)^eM-\[^yl + ^yl}) (15) 
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These hermitian operators can be used to calculate arbitrary expectation values in the 
subsystems. As to the entanglement, assume that |\&) has the Schmidt form (4). Then 

p = i*x*i = E ^K'\K)\^n)(Kmi'\ (19) 

n,n' 

Taking the traces with the |<&") gives n' — n and 

p a = J2\U 2 \K)(K\ , a = 1,2 (20) 

n 

This means that 

• pi and p 2 have the same non-zero eigenvalues 

• these eigenvalues are given by w n = \\ n \ 2 

• their eigenfunctions are the Schmidt functions |$°) 

Therefore the eigenvalue spectrum of the p a gives directly the weights in the Schmidt 
decomposition and a glance at this spectrum shows the basic entanglement features 
of the state, for the chosen bipartition. For this reason, it has also been termed 
"entanglement spectrum" [8]. 

Remarks 

• the p a describe mixed states. An expectation value in subsystem a is given by 

<A a > = J2\^n\ 2 mA a \^} (21) 

n 

• Since the p a are hermitian and have non-negative eigenvalues, one can write 

Pa = \ e- n « (22) 

where Z is a normalization constant and the operator % a has been termed 
"entanglement Hamiltonian" . This form will be encountered permanently in the 
following. 

• The p a should not be confused with e.g. the one-particle density matrices, which 
are simple correlation functions. 

Usually, one starts in a basis where |^) has the form (1). Then 

P=I*X*I= E ^n^, n ,|^)|^)(^,|(^,| (23) 



m,n,m ,n 
r2\ 



and taking the trace with the gives n' — n and 

Pl=EE An,n<J^><^| (24) 



m,m' n 



Thus pi contains the square hermitian matrix AA^ and similarly p 2 contains (A^A)*. 
The form (20) is then obtained by diagonalizing these matrices. This is the general 
approach. 
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Example: Two spins one-half 
A general normalized state is 

|*> = a\+)\+) + b\+)\-) + c|-)|+) + d\-)\-) (25) 

where \a\ 2 + \b\ 2 + |c| 2 + |<i| 2 = 1. The matrix A is then 

*-(:!) <*> 

and one obtains 

\ ca* + bd* cc* + dd* I 
Since the trace is one, the eigenvalues are given by 



wi,2 = \ ± \J\ - det(AAt) (28) 

The state is entangled if Wi >2 ^ 0, 1, i.e. if detA = arf — 6c 7^ 0. This includes the state 
|\I/4) in section 1.3, where a, b, c 7^ and d = 0. 



1.5. Application: DMRG 

The density-matrix renormalization group method (DMRG) is a numerical procedure, 
which was introduced by Steven White in 1992 [9, 10] and makes direct use of the 
Schmidt decomposition and the reduced density matrices. For a review, see [11]. 

Consider a quantum chain, e.g. a spin one-half model, with open ends. Then in 
the simplest variant, the following steps take place, compare fig. 3. 

(0) Start 

Begin with a small system of 5-10 sites. 
Calculate the ground state exactly. 



(1) Schmidt decomposition 
Divide into two halves. 
Calculate the RDM's. 

Diagonalize them and obtain the Schmidt coefficients and Schmidt states. 

(2) Approximation 

Keep only the m Schmidt states with largest weights w n . 
Truncation error: sum of the discarded weights J2 n >m w n- 

(3) Enlargement 

Insert (two) additional sites in the center. 

Form new Hamiltonian in the basis of kept and additional states. 
Calculate ground state. 
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Go back to (1) and repeat. 



(1) 



P 2 



(2) 



m states 




m states 





(3) 



m states 




m states 


• • 



Figure 3. Steps in the (infinite-size) DMRG algorithm. 



For this procedure, the form of the Schmidt spectra is crucial. To have a good 
performance, a rapid drop of the w n is necessary such that only a small number of 
Schmidt states has to be kept. In terms of entanglement, the state must be weakly 
entangled. This is satisfied for non-critical chains. For an Ising model in a transverse 
field, such a small number as 16 Schmidt states gives already a fantastic accuracy for 
the ground-state energy. It is therefore important to understand the features of RDM 
spectra and this leads directly to the study of solvable cases, which is the topic of these 
lectures. 



1.6. Entanglement entropy 

The full RDM spectra give the clearest impression of the entanglement in a bipartite 
system. But it is also desirable to have a simple measure through one number. Since the 
eigenvalues of the RDM's can be viewed as probabilities, one can take the usual entropy, 
as used in probability theory, to characterize the w n . This gives the (von Neumannn) 
entanglement entropy 

S a = -tr(p a In p a ) = w n In w n , (29) 

n 

which is the common entanglement measure for bipartitions. It has the following 
properties 

• Si — S2 = S since the spectra are equal. One can talk of the entanglement entropy. 

• S = for product states. 

• S is maximal if all w n are equal. 

If w n = l/M for n = 1, 2, . . . , M then S = InM. 

The last property leads to a simple interpretation of S. Write 

S = In M eff (30) 
Then e s is an effective number of states in the Schmidt decomposition. 
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A related measure is the Renyi entropy 

S n = -^—\ntr(p n a ) (31) 
1 — n 

where n can also be non-integer. S n has similar properties as S and the same extremal 
values S = and S = InM. For n — > 1 write 

S n = t~~ lntr[p a exp((n - 1) \np a )} (32) 

1 — 71 

and expand the exponential function to obtain Si = S. The Renyi entropy is somewhat 
simpler to calculate, since it contains only a power of p a . The important point is that 
both entropies measure a mutual connection and will, in general, not be proportional 
to the size of a subsystem as usual thermodynamic entropies are. 



1.7. Historical note 

Erhard Schmidt (1876-1959) obtained his PhD in 1905 with Hilbert in Gottingen and 
was professor at the Berlin university 1917-1950. He is most widely known by the 
orthogonalization procedure bearing his name. The work linking him to the quantum 
problems discussed here, appeared in 1907 in the prestigeous journal "Mathematische 
Annalen" [4]. It was based on his thesis and dealt with coupled integral equations with 
a non-symmetric kernel K(s,t). 

In abstract notation, and changing his parameter A to 1/A, the equations were 

Krj} = \<j>, K'<j> = Xip (33) 

He deduced a spectral representation for K 

K(s,t) = J2 X nM*)Mt) (34) 

n 

where <p n and ip n are the eigenfunctions of the symmetric kernels KK' and K'K with 
common eigenvalue 

KK'4> n = X 2 J n , K'Kil> n = \ 2 J> n (35) 

One sees that the kernel K(s,t) corresponds to the total wave function, which for two 
degrees of freedom is \&(a;i, a^)- Moreover, one sees that he already worked with the 
quantities which in the present context are called reduced density matrices. And finally, 
he discussed best approximations for the kernel based on keeping the terms with largest 
weights, which is the same recipe as used in the DMRG. 

The representation of a wave function ^(xi, X2) in this way was discussed in a paper 
by Schrodinger in 1935 [5]. At that time, unsymmetric kernels were already well-known 
in mathematics, so he referred not to Schmidt but to the textbook by Courant and 
Hilbert. The specialists will notice that the equations (33), (35) with K = A — B and 
K' = A + B are just the ones appearing in the famous paper by Lieb, Schultz and Mattis 
(1961) [12] where they diagonalize a quadratic form in fermions. 
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2. Free-particle models 
2.1. Solvable cases 

Before we start the discussion of the free-particle models, which will be the focus of 
these lectures, let me list the quantum states for which one can obtain explicit results 
for bipartite RDM's and thus for the entanglement 

• Ground states of free-fermion or free-boson systems 

• Ground states of certain integrable models 

• Ground states of conformally invariant models 

• Ground states which have matrix-product form or other simple structures 
An example of the last case was the ferromagnetic ground state in section 1.3. 



2.2. Free particles, general result 

Consider models where the Hamiltonian is a quadratic form in fermionic or bosonic 
operators and defined on a lattice. Two standard examples are 

• Fermionic hopping models with conserved particle number 



Coupled oscillators with eigenfrequency uo 

1 d 2 1 



'2 OA + 



(36) 



(37) 



m,n 




For such free-particle models, the reduced density matrices for the ground state can be 
written 



(38) 



Here L is the number of sites in subsystem a and the operators //, f\ are fermionic or 
bosonic creation and annihilation operators for single-particle states with eigenvalues 
Ei. The /'s are related to the original operators in the subsystem by a canonical 
transformation. The constant Z ensures the correct normalization tr(p a ) = 1. 
Note the following features 

• p a looks thermodynamic. 

• the "entanglement Hamiltonian" % a is of the same type as H. 

We will see later that % a is not the Hamiltonian of the subsystem. Therefore (38) is not 
a true Boltzmann formula. Nevertheless, the entanglement problem has been reduced 
to that of a certain Hamiltonian and its thermodynamic properties. But first we want 
to derive the result. 
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2.3. Method 1 - Direct approach 

The direct method to obtain p a is to integrate over the degrees of freedom outside the 
subsystem a. We illustrate it for the example of two oscillators discussed already in 
section 1.3 [13]. The ground state was 

l*o> = (^) 1/4 exp(-iM 1 2 + ^ 2 2 ]) (39) 
One goes through the following steps. 
Step I 

• Write in terms of x\ and x 2 

• Form |^o>(^o|, i.e. ^ (x u x 2 )^ (x[, x' 2 ) 

• Set x' 2 = X2 and integrate over x 2 

• Use (xi + x[) 2 = 2(xj + X?) - {xi - x[) 2 

• Result 

pi(xi,a;i) = C exp (~(a - h)x\) exp {~{x x - x[) 2 ) exp (~(a - b)x' 2 ) (40) 
where a = + u>2)/2 and b = — U2) 2 /2(u>i + u 2 ). 

Due to the derivation, pi has the form of an integral operator. To obtain its 
eigenfunctions and eigenvalues, one would have to solve an integral equation. 

Step II 

• Determine the differential operator for which (40) is the matrix element. 

• Observe that 

exp - x\) 2 ) = 2 {X1 \ exp (1 ^ ) \x[) (41) 

Proof: Express the operator on the right in terms of its eigenfunctions ipk{x) = 
(2"7r) -1 / 2 exp(ikx) and integrate over k 

• Introduce new coordinates y 2 = bx\/2 and the frequency w 2 /4 = (a — b)/b 

• Result 

1 Id 2 1 

Pl = K exp (--wV) exp (- — ) exp (--w 2 y 2 ) (42) 

If one could simply pull the exponentials together, one would have the Hamiltonian of 
a harmonic oscillator in the exponent. However, the exponentials do not commute. 

Step III 

• Write in terms of boson operators a, a* where a = ^co/2(y + 1/ud/dy) 

• Set up equations of motion for Heisenberg operators of a, cc* formed with pi 

• Find Bogoliubov transformation to new boson operators /3, ft 

^cMa + sMat, /3 f = sh^a + ch^a t (43) 

such that becomes a single exponential. This amounts to another stretching of 
the coordinate y — > z. 
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• Result 

Pi — K exp (-e /3 f /3) (44) 

This is the form announced above. The Hamiltonian 7ii in the exponent describes an 
oscillator with frequency e where 

(45) 



, /a ! 

coth(-) = 



'2 y V a-& 2 LV ^2 V wi. 
and its eigenfunctions are those quoted in the Schmidt decomposition (16) when 
expressed in terms of x±. 

This derivation can be generalized to any number of oscillators in a larger system, 
which proves the general statement for this case. However, one sees that the calculation 
involves a number of steps and is already somewhat tedious for the simple case treated 
above. It is therefore fortunate that another much simpler approach exists which we 
will discuss for fermions [14, 15]. 

2.4- Method 2 - Correlation functions 

Consider a system of free fermions hopping between lattice sites with Hamiltonian (36). 
The ground state is a Slater determinant describing the filled Fermi sea. In such a state, 
all many-particle correlation functions factorize into products of one-particle functions. 
For example, 

(4n c l c kCi) = (4Q)(4c fc ) - (c4c fe )(4q) (46) 

If all sites are in the same subsystem, a calculation using the reduced density matrix must 
give the same result. But this is guaranteed by Wick's theorem if p a is the exponential 
of a free-fermion operator 

L 

p a = K exp (- KA c j) ( 47 ) 

where % and j are sites in the subsystem. Thus p a is of the type given in (38). 
The hopping matrix hij is then determined such that it gives the correct one-particle 
correlation functions C it j = (c\cj). This is done in the common diagonal representation 
of both matrices. 

If 4>i(i) are the eigenfunctions of C in the subsystem with eigenvalues Q, the 
transformation 

i 

makes the one-particle function diagonal in the new operators f\ 

(///*'} =G<V (49) 
To obtain this by taking the trace with p a , the operator T-L a must have the diagonal 
form given in (38) with the two eigenvalues related by 

£l = ]n(lzA) or $ = (50) 
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Features 

• Derivation is very short and clear 

• Valid for any Slater determinant 

• Gaussian nature of the problem, only simplest correlator enters 

• Similar for bosonic case 



2.5. Example 

Ring with iV sites and nearest-neighbour hopping. The single-particle states are plane 
waves and H is diagonalized by putting 

Cn = — U eX P (^ n ) C <? ( 51 ) 



N 

q 



In the ground state, the states are filled up to qp and the correlation function is 

C m , n = T7 5^ exp H?( m - n)) {c\c q ) (52) 



N 

q 



/IF A n 
_1 e -"i(™-n) ^ N^oo (53) 
■IF 271 

= sm(q F (m-n)) ^ 
ir(m — n) 

Due to the translation invariance, it depends only on the difference m — n. For half 
filling q F = ir/2. Note the oscillation and the power-law decay of the correlations corre- 
sponding to a critical system. Mathematically, it is a sort of Hilbert matrix. 



Choose a segment of L consecutive sites as subsystem, diagonalize the matrix 
numerically and order the eigenvalues according to their magnitude. This gives the 
following fig. 4. 




Figure 4. Density-matrix spectra for a segment of L sites in an infinite hopping 
model. Left: Single-particle eigenvalues si. Right: Total eigenvalues w n . From [16]. 
Copyright Springer- Verlag, reprinted with permission. 
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Features 

• Dispersion of S\ roughly linear with curvature 

• Values of order 1 and larger 

• Curves flatter for larger L 

• Rapid initial decrease of the w n 

• Entanglement small, but increasing with L 

2.6. Characteristics of the problem 

(a) Single-particle eigenfunctions 

For the low- lying the eigenfunctions are localized near the boundaries. 
This is shown in fig. 5 for a non-critical and a critical hopping chain. 




20 40 60 80 100 20 40 60 80 



j j 

Figure 5. Lowest lying single-particle eigenstates in a dimerized (6 = 0.1, left) and a 
homogeneous (S = 0, right) hopping model for a segment of L = 100 sites. From [16]. 
Copyright Springer- Verlag, reprinted with permission. 



Consequences 

• Double degeneracy of low E\ for segments in non-critical chains 

• Slower decay of the resulting w n 

In fig. 6 this is illustrated for a half-chain of coupled oscillators. Shown are the results 
both for the open chain (fig. 1 left) where one has one boundary and for the ring (fig. 
1 centre) where the subsystem is a segment with two boundaries. 

The slower decay of the w n leads to a poorer performance of the DMRG for rings 
and explains why the method is normally used in the open-chain geometry. In two 
dimensions, whole bands of E\ arise which are associated with the boundary between 
the subsystems, see section 4.4. In the entanglement entropy this leads to the so-called 
"area law" . 




(b) Entanglement Hamiltonian T-L a 

In general, this operator is different from the Hamiltonian of the subsystem. 

This is shown in fig. 7 for a segment in a hopping chain. The hopping matrix hij 
in (47) was calculated, using the common eigenfunctions 0/ of C and h, via 

h ij = ^2<f>i(i)£i<f>iU) ( 55 ) 
i 




2 4 6 8 10 12 14 16 



Figure 7. Matrix elements in "H Q for a hopping model. First, third and fifth neighbour 
hopping in a segment of L = 16 sites. From [1]. Copyright IOP Publishing, reprinted 
with permission. 



The dominant elements are those for nearest-neighbour hopping and vary roughly 
parabolically, whereas in the chain they are constant. For a half-chain one finds half a 
parabola. 
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(c) Spectrum of C 

In a large subsystem, most of the eigenvalues Q he (exponentially) close to and 

to I. 

This is illustrated in fig. 8 for a segment in a hopping model 

i ...... 

0.8 
0.6 - 

tl 

0.4 
0.2 - 

o 1 ^ 

2 4 6 8 10 12 14 16 18 20 

1 

Figure 8. Eigenvalues of the correlation matrix for a segment of L = 20 sites in a 
hoppping chain. 

It can be understood from (52) as follows 

• In the total system, the eigenvalues are {c q c q ) = 0, 1 

• Restricting C to the subsystem changes the spectrum 

• Low-lying states localized near the boundary appear, compare (a) 

• But bulk states remain 

In numerics, this leads to the following difficulty. The closeness of Q to or 1 soon 
exceeds the usual double-precision accuracy. The Si can then no longer be determined 
reliably, unless one works with special techniques. Therefore the values of the Si in most 
of the figures shown here do not exceed 20-30. However, for the entanglement this does 
not matter, since large Si give negligible contributions. 

A special role is also played by eigenvalues Q = 1/2 corresponding to S\ = 0. Such 
an eigenvalue causes a two-fold degeneracy of all w n and is therefore seen in the RDM 
spectrum. These zero modes have found much interest recently because they may reflect 
a symmetry of the real Hamiltonian with boundaries. 

2.7. Schmidt form for fermions 

The correlation function approach gives the Schmidt spectra in a very easy way. But 
it is also instructive to derive the Schmidt decomposition directly. This is done in the 
following way [18]. 

• Consider a system with N particles. Divide the occupied single-particle states tp q (n) 
into the components 

I iPq{n) : n e 2 
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These are neither orthogonal nor normalized in their subsystems. 

• Find new states Xi( n ) sucn that their components xl( n ) an d xf ( n ) are orthogonal 
in their subsystems. This is done by diagonalizing the overlap matrices 

Ks' = a = 1,2 (57) 

Their eigenvalues are Q and 1—0 and the new functions have the norms 

{Xi\Xi) = Ch (XilXi> = l-0 (58) 

• Form normalized states via 

4>\ = ^xl = 7yf=^ ( 59 ) 

• Define Fermi operators a a j for the <f>f. Then 

N 

l*> = II + V^O^i] |0> (60) 

1=1 

where |0) is the vacuum. This gives the Schmidt decomposition if one multiplies 
out the product. 

Comments 

• Instead of the L x L correlation matrix C, the N x N overlap matrix M appears 

• However, the non-trivial eigenvalues Q are the same 

• A particle in state Xi is found with probability Q in part 1 and with probability 
1 — Q in part 2 

• If Q — the particle is found only in subsystem 2. This has to happen, if subsystem 
1 cannot accomodate all the iV particles. 

• The approach can be applied to continuous systems where ip q {n) — > ip q (x) 

The approach shows that the single-particle eigenvalues Q in one subsystem are associ- 
ated with the eigenvalues 1 — Q in the other. The two lead to ±ei and give the same 
w n -spectrum, as it should be. 

Example [19] 

iV free fermions on a ring of length L, subsystem segment (—£/2,£/2). 
Single-particle wavefunctions 

1 2% 

ip q (x) = —= exp (iqx), q = — n, n — 0, ±1, ±2, . . . (61) 
v L L 

Overlap matrix in subsystem 

^ = WW) (62) 
1 f e/2 

— — dx exp (—i(q — q')x) (63) 
L J-e/2 

2 — sin((q-q')£/2) (64) 



(q - q')L 
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Writing q = 2nm/L, q' = 2nn/L, the matrix becomes 

M^ = Sin((7r y L)(m ~ n)) (65) 

This is the correlation matrix result (54) with the substitution qp — > n£/L. The case 
£ = L/2 corresponds to half filling, and one can take over the lattice results for the (i. 
Choosing a different segment of the same length changes M but not the eigenvalues. 



2.8. Some additional details 

• In the correlation function approach, the eigenvalue equation can also be written 
in the form 

(l-2C)^ = tanh(|)0,. (66) 

• If the expectation values F it j = (cjcp and F*j = (cjCj) are non-zero, they have to 
be included in the considerations. Then for real F the equation becomes 

(2C - 1 - 2F)(2C - 1 + 2F) fa = tanh 2 (^) fa. (67) 

• Instead of working with the usual fermions, one can use Majorana fermions defined 
by 

a 2 n-l = (Cn + 4), 02n = i(Cn-4) ( 68 ) 

and form the 2L x 2L correlation matrix {a m a n ) in the subsystem. It has eigenvalues 
1 ± ?tanh(^/2). This is usually done if the "anomalous" correlation functions F i: j 
exist. 

• For coupled oscillators, the correlation functions of position variables and of 
momenta, X it j = (xiXj) and Pjj = (piPj), take the place of the Majorana variables. 
The single-particle eigenvalues then follow from 

2P2X^ = coth 2 (|)0,. (69) 

For the two coupled oscillators treated in section 2.3, one has (x\) = {I / uoi + l / uo 2 ) / A 
and (p\) = (ooi + u; 2 )/4 which gives again (45). 
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3. Integrable models 

In one dimension one can exploit the relations between quantum spin chains and 
two-dimensional classical models. For non-critical integrable models, this allows to 
determine the RDM's and their spectra analytically for large systems divided in the 
middle. 

3.1. Transverse Ising model 

We will discuss the the approach for the Ising model in a transverse field (TI model) 
with Hamiltonian 

# = -E<- A E«+i' ( 7 °) 

n n 

The transverse field has been set to h — 1. The ground state is non-degenerate for A < 1 
and asymptotically degenerate with long-range order for A > 1. If rewritten in terms of 
Fermi operators, H becomes a quadratic form 

H = - £(2 4 Cn - 1) - A - c n )(cl +1 + c n+1 ). (71) 

n n 

Therefore, according to section 2 

1 L 
Pa=z e ~ Ha > U a = Y,tif}fi (72) 

i=i 

and the E\ could be calculated numerically using the correlation functions. The present 
approach will give them analytically. 

3.2. Relation to a 2D partition function 

The TI model has the following features 

• H commutes (up to boundary terms) with a particular (diagonal) transfer matrix 
T of an isotropic 2D Ising model on a square lattice 

• Its ground state is the eigenstate of T with maximal eigenvalue 

From the second property, it follows that one can obtain |\1>) from an initial state \^ s ) 
via 

~ lim T n \^ s ) (73) 

n— >oo 

In this way, one has related to the partition function of a two-dimensional semi- 
infinite Ising strip. This is a discrete path-integral representation of It follows 
that 

• p = l^) is given by two such strips 

• p a is obtained by tying the two half-strips together 

In this way, p a is expressed as the partition function of a fully infinite strip with a 
perpendicular cut. This is shown in Fig. 9 on the left. 
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Figure 9. Left: Density matrices for a quantum chain as two-dimensional partition 
functions. Far left: Expression for p. Half left: Expression for p\. The matrices 
arc defined by the variables along the thick lines. Right: Two-dimensional system 
built from four quadrants with corresponding corner transfer matrices A, B, C, D. 
The arrows indicate the direction of transfer. From [16]. Copyright Springer- Verlag, 
reprinted with permission. 



3.3. Some transfer matrix formulae 

Before we discuss the evaluation of this particular partition function, we list a few rela- 
tions for conventional Ising transfer matrices. 

(a) One dimension 

Consider the Ising chain with Hamiltonian 

H = -Jj3<7„<7 n+ i, (74) 

n 

where a n = ±1. To calculate a partition function, one needs (K = (3J) 
exp (—/3H) = exp (Kaia 2 ) exp (Ka 2 a 3 ) exp (K0-3G4) . . . 

= T((Ti, <j 2 ) T(<7 2 , 0-3) T(a 3 , cr 4 ) . . . (75) 
Each T contains the Boltzmann factor for one bond and is a 2 x 2 matrix 



e K e~ K 
e~ K e K 



(76) 



Summing over all a n = ±1 multiplies the matrices together and gives, for a ring of N 
sites, the partition function Z = trT N . In operator form, T can be written 

T = C exp (K*a x ) (77) 

with the so-called dual coupling K* defined by sinh2.fr* = 1/ sinh2i^. It is large, if K 
is small and vice versa. 

(b) Two dimensions 

In two dimensions, one can build up a lattice row by row. The transfer matrix then 
contains the Boltzmann factors for the vertical and horizontal bonds in one row. This 
is shown in fig. 10 (a) by the thick lines. 
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Figure 10. Geometry for two types of transfer matrices, (a) Row-to-row transfer 
matrix (b) Corner transfer matrix. 

For N sites in a row, T is now & 2 N x 2 N matrix and in operator form given by 

T = T\T 2 = C N exp (K{ £ <) exp (K 2 £ (78) 

n n 

where K\ and K 2 are the vertical and horizontal couplings, respectively. 
Features 

• Terms like in transverse Ising model 

• However, the exponentials do not commute 

• Exception: K*,K 2 1, strong vertical and weak horizontal bonds. Then one can 
combine the exponentials. This is called the "Hamiltonian limit". 

3.4- Corner transfer matrices 

To calculate the partition function needed for p a , a kind of "circular" transfer matrix 
would be appropriate. This is indicated in fig. 9 on the right. Then p a would be given 
by 

p a ~ ABCD (79) 

It so happens that such quantities were introduced by Baxter in 1976, see [20]. It turned 
out that for integrable models they have fascinating and simple properties which make 
them a powerful tool for calulating order parameters. 
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(a) Structure 

In fig. 10 (b) a quadrant of a square lattice model is shown. The CTM contains 
all Boltzmann factors indicated by thick lines. The internal variables are summed. In 
operator form, leaving out the prefactor C 

• Horizontal bonds 

exp (K 2 o-qO-1) (1), exp (K 2 a\a z 2 ) (3), exp (K 2 o z 2 al) (5), . . . 

• Vertical bonds 

exp(fl>f) (2), exp of) (4), exp (K* of) (6),... 

• All matrices to be multiplied in correct order from bottom to top 
Hamiltonian limit 



with 



A = e~ n °™ (80) 
Uctm = K{ 2n < + K 2 £(2n - 1) « +1 (81) 



n>l n>l 

Features 

• Inhomogeneous TI Hamiltonian 

• Fields and couplings increase linearly 

• Eigenvalues equidistant for L — > oo 

F _ f (21 ~ l)e , K* < K 2 

£l ~\ 2le , K*>K 2 (82) 

• To be seen directly in the limiting cases 

• Otherwise result of a fermionic calculation 

(b) General case 

So far only the Hamiltonian limit has been considered. The structure of Uctm is 
then a consequence of the wedge-like geometry. However, for determining p a via (79), 
this is not enough, since in the next quadrant the anisotropy is the other way around. 

Amazingly, however, the following holds asymptotically 

• The eigenvalue spectrum of %ctm has the form (82) for arbitrary couplings 

• In the product ABCD, the parameter e which gives the level spacing is 

e = irl(k')/l(k), (83) 

where k with < k < 1 is either given by k = sinh 2K\ sinh 2K 2 or by 
k — 1/ sinh2i^! sinh2i^ 2) whichever is smaller than 1. I(k) is the complete elliptic 



integral of the first kind and k' = yl — k 2 . 
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Figure 11. Level spacing as a function of the parameter k. From [1]. Copyright IOP 
Publishing, reprinted with permission. 




The parameter e diverges for k — > and vanishes for A; — > 1, in both cases 
logarithmically. It is shown in fig. 11. 

The derivation uses the integrability of the model, which is contained in the so-called 
star-triangle equations, and a proper elliptic parametrization of the couplings. This 
leads to two parameters, the k appearing above which is connected with the tempera- 
ture, and another parameter u which measures the anisotropy, but does not enter the 
product ABCD. A brief account can be found in the Les Houches lectures of Cardy 
1988 [21]. 

(c) Application to RDM 

The CTM discussed so far can be used for calculating the spontaneous 
magnetization as expectation value of the central spin. This is sketched in the 
supplement. However, this central spin is an obstacle for the RDM application, because 
it is common to all four CTM's and prevents the division of the system into two parts. 
To calculate the partition function for p a , one uses the modified CTM shown in fig. 12. 
This amounts to an interchange of the coefficients 2n and 2n — 1 in (81), (82). 




Figure 12. Corner transfer matrices without central spin for calculating the RDM. 
Left: single matrix. Right: arrangement of four such matrices giving p a . 
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Summing up, the result for the single-particle eigenvalues in % a is 

ei = 



(2Z + l)e 
2le 



disordered region 
ordered region 



(84) 



where I = 0, 1, 2, . . . and e is given by (83). In terms of the TI model, the parameter k 
is 



k 



A 
1/A 



A < 1 
A > 1 



(85) 



3.5. Spectra and entanglement 

In Fig. 13, spectra are shown for a finite open TI chain with N = 20 sites, divided in 
the middle. Thus the subsystem has L — 10 sites and there are 10 eigenvalues £/. The 
example displays both the infinite-size properties and the modifications by the finite 
size. 

Features 

• Linear behaviour of e\ as predicted 

• Deviations at upper end closer to the critical point A = 1 

• At A = 1 shape as for hopping model 

• w n decrease extremely rapidly for small A (note the scale) 

• uv-decay slower near criticality, but still impressive 




9 10 




Figure 13. Density-matrix spectra for one-half of a transverse Ising chain with TV = 20 
sites in its ground state. Left: All ten single-particle eigenvalues Right: The largest 
total eigenvalues w n . From Chung [17]. 

This means that the ground state is weakly entangled. A Schmidt decomposition 
can be truncated safely after about 10 terms and this is the explanation for the fantastic 
performance of the DMRG in this case [22]. Note that altogether there are 2 10 = 1024 
w n already in this small system ! 
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Behaviour of the w. 



n 



• Plateaus in w n for strictly equidistant levels 

• Behaviour for large n from number of partitions 

w n ~ exp[— a(lnn) 2 ] 



(86) 



where a = e 6/n 2 . 

The case of a segment cannot be treated by the CTM method, but one can simply 
include the degeneracy seen numerically (section 2.6) into the CTM results. Segments 
in free-particle models can be treated by a different method which, however, is more 
technical and less physical than the CTM approach [23]. 

3.6. Other systems 

The CTM approach works also for a number of other quantum chains, namely 
• The XY spin chain with Hamiltonian 



This generalization of the TI chain also corresponds to a free-fermion problem. 
2D problem: Ising model on a triangular lattice 

• The XXZ and XYZ Heisenberg spin chains which contain fermion interactions. 
2D problem: Eight-vertex model 

• The oscillator chain with nearest-neighbour coupling. 
2D problem: Gaussian model 

It turns out that the CTM spectrum has the form (84) for all these models, even if 
they contain interactions. Thus one has a universality in these problems which makes 
the entanglement properties of all the fermionic systems identical. Only the parameter 
k is related differently to the system parameters in each case. For the oscillator chain, 
for example, it is given by k/(l — k) 2 = K/oo 2 if K is the nearest-neighbour coupling. 
This chain is the bosonic analogue of the TI chain, but it has no ordered phase. In spite 
of the different statistics, the w n spectra are similar and the asymptotic law (86) holds 
with a smaller a. 

The bosonic formula can also be used to treat exactly a two-dimensional lattice of 
coupled oscillators which is divided in the middle by a straight line. This is because by 
making a Fourier transformation parallel to the interface, the problem separates into 
uncoupled chains. 

3.7. Supplement: Onsager formula 

With the CTM spectra for the 2D Ising model, the famous Onsager formula for the 
spontaneous magnetization can be derived in a few lines. Working in the geometry of 




n 



n 



(87) 
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fig. 10 (b) and fixing the outer spins as indicated, the expectation value of the central 
spin has the form 

W> = (88) 

where Z + and Z_ are the partition functions with o"o parallel and antiparallel to the 
boundary spins, respectively. In terms of the CTM's, this becomes a quotient of traces 

tr(a^iABCD) 

{ao) = tr(ABCD) m 
In the fermionic representation, the operator o-§o z L can be expressed in terms of the 
operators which diagonalize Hctm as exp (in J2i fi fi)- The trace can then be performed 
for each / separately and the exponential factor leads to a minus sign in the numerator. 
Thus 

Since one has to consider the ordered region, one has to choose e t = (21 — l)e in (82). 
With q = e~ £ the product then is 

i=i H 

Due to its definition, q is an elliptic no me and the infinite product in (91) has a simple 
relation to the elliptic moduli k and k' which appear in e. This gives 

(a ) = (fc') 1/4 =(l-fc 2 ) 1 / 8 (92) 

which is Onsager's formula. The parameter k is here k = 1/ sinh2i^! sinh2i ; r 2 . It is 
interesting to note that also Yang in his 1952 proof of Onsager's result [24] derived an 
infinite product equivalent to (91), although his approach was quite different. In the 
CTM formalism, it appears in a natural way, and also the order parameters for more 
complicated models take such product forms, see Baxter's book. 
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4. Entanglement entropies 

We have seen already some RDM spectra, which contain the full entanglement 
information. In this section we want to see how their properties translate into the 
entanglement entropy. Entanglement entropies are the standard quantities considered 
in this area and have been the topic of a large number of studies. 

4-1. General 

Due to the form of the p a , one has the same expressions for the von Neumann entropy 
as in thermodynamics. Thus, F = U — TS with T — 1, or S — —F + U, and the 
free-particle character of T-L a gives, as in statistical physics 



where the upper (lower) sign refers to fermions(bosons). From this formula, one can 
immediately see some general properties 

• Largest contributions come from small E\ 

• Therefore entropy particularly large in critical systems 

• Maximum value for fermions L In 2 if all E\ = 

• If all ei are m-fold degenerate, S has m times the value without the degeneracy 

The last property is an additivity which appears e.g. for uncoupled chains or for 
two independent interfaces. As to the magnitude, an eigenvalue ei ~ 1 also gives a 
contribution of order 1 to S and the sums converge rapidly for larger S\. 

4-2. Example: TI chain 

With the spectra found in section 3, it is easy to calculate S for the infinite transverse 
Ising chain. The result of a numerical evaluation is shown in fig. 14 

One notes the following features 

• S vanishes for A — > 0. 

Formally: All e t diverge. Physically: |^) becomes product state. 

• S goes to In 2 for A — > 00. 

Formally: All Si except one diverge, e is zero. Physically: |\&) is superposition of 
the two product states | + + + . . .) and | . . .}. 

• S diverges at the critical point A = 1. Formally: e — > 0, slope of the dispersion 
curve goes to zero. Physically: State becomes more and more entangled as the 
correlation length increases. 

Due to the equidistant single-particle levels, one can even calculate S in closed 
form. In the disordered region, one finds with k = A 




(93) 




(94) 
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Figure 14. Entanglement entropy between the two halves of an infinite TI chain as a 
function of A. From Calabrese and Cardy [25]. Copyright IOP Publishing, reprinted 
with permission. 



A similar expression with an additional contribution of In 2 coming from the eigenvalue 

Eq = holds in the ordered region. 

From this, one can extract the behaviour near k = 1 

and since the correlation length is given by £ ~ 1/(1 — k), this can be written 

S=^rn£ (96) 

which shows a logarithmic critical behaviour. The effective number of states in the 
Schmidt decomposition, however, has normal power-law behaviour 

M eS ~ e /U (97) 

In this sense, the coefficient of the logarithm is a critical exponent. 
The Renyi entropies are 

1-n^ (l + e- £ O rt 
and lead to more complicated closed expressions, but the critical behaviour is analogous 

Sn = ^(l + -)hi£ (99) 
24 n 

An unusual structure is seen if one looks at the next (subleading) terms in the expansion. 
One finds that they are of the form £ _fc / n with k = 1, 2, 3 . . ., i.e. the powers depend on 
the Renyi index n which determines the number of windings in the path integral for p™ 
[26, 27]. The same phenomenon is encountered for order parameters on such Riemann 
manifolds. 
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4-3. Critical chains 

At a critical point, one has to work with finite subsystems. The spectra for a hopping 
model have already been shown in section 2.5, and a marked size dependence was noted. 
The dispersion curves of the e t became flatter with increasing L. This gives an increase 
of S. From (96) one can already guess that £ will be replaced by the length of the 
subsystem, and this is in fact the case. The asymptotic formula is 

S = u^lnL + k (100) 

Features 

• v — 1,2 number of contact points between subsystem and the rest 

• k non-universal constant (subleading term) 

• c central charge, from conformal considerations, c = 1/2 for TI model, c = 1 for 
hopping model 

This result can be understood for the hopping model as follows. The Si curves for 
small systems are not linear, but show curvature. However, for large L, more precisely 
for large InL, one can use a continuum approximation to the eigenvalue equation to 
derive the formula, for a segment in a chain, 

£l = ± 2lnT (2/-1) ' / = 1 ' 2 ' 3 --- ( 101 ) 

Using this in (93) and changing the sums into integrals gives 

2\nL \ f°° . „ f°° J e 

S = — / de In 1 + exp(-e)) + / de — - — 

[Jo Jo exp(e) + l 

and since both integrals equal 7r 2 /12 one finds 

S=\ InL (103) 

In numerical calculations, this logarithmic law can be seen already in relatively small 
systems, where (101) does not yet hold, but one has approximately InL — > InL + 2.5 
for the first eigenvalues. 

The expression for the Renyi entropy follows in the same way by going over to 
integrals in (98) and gives for a segment 

S n = hl + -)\nL (104) 
o n 

4-4- Higher dimensions 

As mentioned in section 2.6, one finds bands of E\ in two dimensions. This is illustrated 
in fig. 15 for a 10 x 10 square lattice of coupled oscillators, divided into two halves. 
The vertical coupling was varied and one can see how the plateaus with 10 levels (for 
10 uncoupled chains) develop into bands. The states in a band can be indexed by a 



(102) 
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Figure 15. Single-particle eigenvalues for one-half of a 10 x 10 system of coupled 
oscillators with uj = k x — 1 and different couplings k y . From [28]. Copyright APS, 
reprinted with permission. 



vertical momentum q y = q. 

For the entropy, this has the following consequences 

• Without coupling: each chain gives the same contribution s to the total 
entanglement entropy. Thus for M chains one has S = M s. 

• With coupling: one has to add up the contributions s(q) for each value of q. For 
large M 

S = J2s(q)^M J*^s(q) (105) 

• Therefore S proportional to the length of the interface between the subsystems. 

• In three dimensions: area of the interface 

• Also for other geometries 

This is the so-called area law for the entanglement entropy. For fermionic critical 
systems, however, one has logarithmic corrections. For a system with typical size L in 
d dimensions, one finds 

S ~ L d_1 In L (106) 

if the the state corresponds to a finite Fermi surface. This can be proved exactly by 
putting bounds on S [32, 33], see section 4.7. 

4-5. Entanglement across a defect 

Since the entanglement is a kind of boundary phenomenon, one expects that in will 
be changed by a modification of the interface between the subsystems. This has been 
investigated for hopping chains and critical TI chains with a modified bond, as shown 
in fig. 16. 
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Figure 16. Transverse Ising chain with a bond defect. 



Limiting cases 

• Chain cut by defect, t — 0: no entanglement, S = 

• Chain homogeneous, t — 1: logarithmic law (100), S ~ InL 

What happens in between ? Numerical results for the E\ are shown in fig. 17. 
Features 

• Development of a gap at the lower end of the spectrum 

• Upward shift of the whole dispersion curve as t goes to zero 

• Therefore decrease of S for fixed L 

• Logarithmic law for S remains valid 

• But c — > c e g(t). 

The variation of c e fr with t can be determined numerically, but it turns out that it 
can also be calculated analytically [29]. Since it is an exercise in going to two dimensions 
and using partition functions as in section 3, it is presented here briefly. Because one is 
at the critical point, one can use conformal mappings. The scheme is shown in fig. 18. 

In the end, one obtains an expression for the Si with a gap which one can insert into 
the continuum formula (102). The integrals lead to dilogarithms in terms of a parameter 
s — 2/(t + l/t) which is the transmission amplitude through the defect, i.e. s 2 is the 
transmission coefficient. The formula is somewhat long, so it is more instructive to show 
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Figure 17. Single-particle eigenvalues £; as a function of the defect strength for 
TI chains with 2L = 300 sites. From [29]. Copyright Wiley- VCH, reprinted with 
permisson. 
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Figure 18. Representation of p a for a chain with a defect in the centre by two- 
dimensional partition functions. Left: Original representation. Centre: Simplified 
annular geometry. Right: Strip geometry obtained via the mapping w = Inz. The 
defect line is always shown dashed. From [29]. Copyright Wiley- VCH, reprinted with 
permisson. 
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Figure 19. Effective central charge c e g(t) for a TI chain as a function of the defect 
strength t (left curve). From [29]. Copyright Wiley- VCH, reprinted with permisson. 

the result graphically, see fig. 19. For the Renyi entropy S 2 , by the way, one finds a 
very simple result, namely 

c c ff,2 = ^■arcsin 2 (s/v / 2) (107) 

7T 

The continuous variation of the coefficient might seem natural, but it is connected 
with the free-fermion nature of the TI and the hopping chain. The defect is then a 
"marginal" perturbation which changes also the local magnetic exponent continuously. 
Things are different for a defect in an XXZ chain, which is a Fermi system with 
interactions. Then a defect either leads to c e g = if the interaction is repulsive, or 
is irrelevant, i.e. c c g- = 1, if the interaction is attractive. This is in analogy to the 
transmission properties in this case. 
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4-6. Inhomogeneous systems 

The entanglement can decrease or increase if one makes a system inhomogeneous. This 
is illustrated here with two simple but instructive examples. 

(a) Hopping chain in a field [30] 

Consider an open chain of 2L sites with Hamiltonian 

^ L-l L 

H = — 2 E (clc n+1 + cl +1 c n ) + h J2 (n-l/2)cic n (108) 

n=-L+l n=-L+l 

This describes the so-called Wannier-Stark problem of electrons in a constant electric 
field. In magnetic language, it is an XX chain with a linearly varying magnetic field in 
the z-direction. Due to the field, the particles accumulate on the left. 

Features 

• Density profile, system full on the left and empty on the right 

• Characteristic length A = 1/h 

• Transition region has width 2A 

• Single-particle wave functions are Bessel functions <frk{n) = 
concentrated near site k. 

• Single-particle energies are equidistant, Uk = h(k — 1/2), Wannier-Stark ladder 

Correlation matrix for a half-filled system for L — > 00 

00 

C mn = ^ Jk+m{tyJk+n{ty (109) 
fc=0 

= 777 r [Jm-lJn — JmJn-l] (HO) 

2(m — n) 

In the limit A — > oo, this reduces to the result (54) for the homogeneous chain. 

The length scale A is seen also in the low eigenvectors of C. They are essentially 
confined to the transition region. 

Numerical results for the entanglement entropy if the system is divided in the middle 
are shown in fig. 20. 

Features 

• Logarithmic up to L A 

• Saturation for L > A, if h ^ 

• Saturation value for large A 

S oa {\) = \\n{2\) (111) 
o 

This is analogous to (96), where the correlation length entered. Interpretation: The 
parts outside the interface region, which are either full or empty, cannot contribute to 
the entanglement. 
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Figure 20. Entanglement entropy for a hopping chain in a linear potential as a 
function of the half-length L. From [30]. Copyright IOP Publishing, reprinted with 
permission. 



(b) Inhomogeneous hopping [31] 

Consider a model with Hamiltonian 

1 

H =~2 22 tn(clc n+1 + ci +1 c n ) (112) 
n=-L+l 

where the hopping amplitudes t n decay rapidly from the center towards the ends of the 
chain, for example like t n = exp(— |n|). In this model, the density in the ground state 
is constant as for a homogeneous chain. However, the state is highly entangled. 



Example: Four sites 



h h 



Figure 21. Four-site chain with corresponding hopping amplitudes. 



For ti <C t , the lowest single-particle states have energies u± = —t and 
u 2 = —tf/to- These states are occupied in the ground state and the corresponding 
eigenvectors are approximately 
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In the first one, sites 2 and 3 are fully entangled, in the second one sites 1 and 4. 
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The total correlation matrix is 



°"5 



/ 



V 



1 

1 1 

1 1 

-10 



-1 \ 






(114) 



Restricting C to the left or right half-chain, one finds Ci = C2 = 1/2, i.e. E\ = e 2 = 0, 
which gives S = 2 In 2. The mechanism persists for larger systems and leads to the 
concentric structure shown in fig. 22. 



Figure 22. Concentric entanglement structure in an inhomogeneous hopping model. 
After Vitagliano et al. [31]. 



4-7. Entropy and fluctuations 

In hopping models, there is a close connection between the entanglement entropy and 
the particle-number fluctuations in the considered subsystem. This allows to put bounds 
on S [32, 33]. 

In terms of the eigenvalues Q of the correlation matrix C, one has 

S = - [Ci In + (1 - Ci) ln(l - 0)] = E s ^ ( 115 ) 
1 1 

The function s(x) defined by minus the bracket in (115) has the properties 

• Symmetry with respect to x — 1/2 

• s(x) = for x = and x — 1 

• Maximum at x — 1/2 with s(l/2) = In 2 

As a result, it can be bounded in < x < 1 by a parabola 

s(x) > 4m2 x(l -x) (116) 
and the equality holds for x — 0, 1/2, 1. This is shown graphically in fig. 23 
It follows that 

S > 4 In 2 E0(1 " 0) = 4 In 2 tr[C(l - C)] (117) 

But the traces can be written as 

tr[C(l-C)] = (iV 2 >-(iV> 2 (118) 

where iV = c\ci is the particle number operator in the subsystem. Therefore the 
particle-number fluctuations give a lower bound on S 

S > 41n2 [(N 2 ) - (N) 2 ] (119) 
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X 

Figure 23. The function s(x) (solid) and its quadratic lower bound (dashed). 



These fluctuations have a direct physical significance and are easier to calculate. 
Application 

• One dimension, large L 

[(N 2 )-(N) 2 ] = ij\nL (120) 

• Two dimensions, large L 

[(N 2 ) - (N) 2 } ~ L InL (121) 

By shifting the parabola x(l — x) upwards, one can also obtain upper bounds. In this 
way, one can prove the behaviour of the entropy in various dimensions without actually 
calculating it. The lower bound in ID gives the prefactor 41n2/7r 2 = 0.28, which is 
rather close to the exact value 1/3. From these considerations, a general formula for the 
prefactor was obtained which involves an integral over the surface of the subsystem in 
real space and the Fermi surface in momentum space [33]. 
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5. Quenches and miscellaneous 

So far we have been concerned with time-independent situations. In this last section, 
we turn to cases where the entanglement changes in time. Moreover, I return once more 
to possible relations between the entanglement Hamiltonian and the real one and finally 
give a short summary. 

5.1. Quenches 

If a quantum state changes in time, this will in general affect the entanglement 
properties. However, the change must be more than a mere phase factor. Thus one 
has to have a time evolution with a Hamiltonian, for which is not an eigenstate. 
The simplest set-up is to make an instantaneous change 

H -> H 1 (122) 

After that 

• The state |\1>) evolves as \${t)) = e~ iHlt \^ ). 

• The total density matrix p evolves. 

• The RDM's p a also evolve. 

If Hi is a free-particle operator, the arguments work as before. If the initial state 
was a Slater determinant, the correlation functions at time t 

(*(*)|4cKq|*(0> = (*o|4W4(t)c fc (t)Q(0l*o) (123) 

factor again, because the Heisenberg operators Ck(t) at time t are then linear 
combinations of the initial ones. Therefore p a {t) has the exponential form (38) but with 
a time- dependent operator H a {t) and the eigenvalues Si{t) follow from the correlation 
matrix at time t 

C iJ (t) = (*o\<i(t)c j (t)\* ). (124) 

Therefore, one only needs to determine the time evolution of the operators Cj(t) in the 
Heisenberg picture. 

Physically, one finds a surprising phenomenon, namely the entanglement increases after 
the quench 

• In global quenches S ~ t 

• In local quenches, S ~ \nt 

We show this explicitly for two examples. 
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5.2. Global quench 

Hopping model 

• Start from fully dimerized, half- filled model, only pairs of sites (2n, 2n + 1) are 
coupled and correlated. 

• Make it homogeneous with dispersion relation oj q = — cos q and let it evolve. 



The time evolution of the Fermi operators then involves Bessel functions 

C i(^) ^ ] ft Jj—m{f)C-m 



and the result for the correlation matrix is 



Cm,n(t) — ^ 



S m ,n + \(Sn,m + l + S n , m -l) + e"**^ J m _ n{2 t) 



(125) 



(126) 



The resulting single-particle spectra are shown on the left of Fig. 30. 
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Figure 24. Global quench in a hopping model, starting with a fully dimerized initial 
state. Left: Time evolution of the single-particle spectrum for a segment of L = 100 
sites. Right: Entanglement entropy with the asymptotic value. From [30]. Copyright 
IOP Publishing, reprinted with permission. 



Features 

• Dispersion linear near zero 

• Slope decreases with time, S increases 

• For times t 3> L/2 approach to a limiting curve, S saturates 

The asymptotic form of the spectrum follows from the first three terms in (126) which 
correspond to a tridiagonal correlation matrix and are the Fourier transform of the 
constants (cjc g ) in the initial state. The eigenvalues for a segment are 



1 7T 

Cj(oo) = 2^ + cos <?0> qi = J^l l i l = l,2...L 



(127) 



and lead to 



Ei(oo) = 2 lntan(g//2). 



(128) 
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The spacing of the qi is proportional to 1/L and gives an extensive entropy S = 
L(21n2 - 1). 

The build-up of an extensive entropy is a typical signature of global quenches. 

It has a simple physical interpretation due to Calabrese and Cardy [34] sketched in fig. 

25. 

• Particle-hole pairs are emitted 

• Create entanglement between the subsystem and remainder 

• Travel with maximum velocity v — 1 

• "Light-cone effect" , S ~ t as long as separation 2t < L 




Figure 25. Creation of entanglement after a global quench by emitted particle-hole 
pairs for the case of a segment in a chain. 

The result is relevant for numerical calculations, because it means that one can follow the 
evolution only for a limited time with DMRG. Beyond that, the state is too entangled 
to be well approximated. 

5.3. Local quench 

Hopping model, set-up shown in fig. 26 

......... 

.............. 

Figure 26. Two variants of a local quench. 

• Initially subsystem (center) decoupled from the rest 

• Add bond(s) to create a homogeneous chain and let system evolve 



The evolution of the Fermi operators is again given by (125), but the initial 
condition is different. The calculation has to be done numerically. In fig. 27 the 
result for 5* is shown. 
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infinite 




Figure 27. Entanglement entropy for the two geometries after the quench for a 
subsystem of length L = 40. From [36]. Copyright IOP Publishing, reprinted with 
permission. 



Features 

• "Entanglement bursts" after the connection 

• Duration t = L (infinite case) and t = 2L (semi-infinite case) 

• For larger times approach to equilibrium (dotted) 

The plateau can be related to a front which starts from the initial defect site and 
travels through the subsystem until it leaves it again. This is seen directly in the lowest 
eigenvector in fig. 28. Using methods of conformal field theory, one can derive analytical 




Figure 28. Front propagation in the lowest single-particle eigenvector for the semi- 
infinite geometry and L = 100. Shown are the times t = 20,60, 120, 160. From [36]. 
Copyright IOP Publishing, reprinted with permission. 



formulae for both cases [35, 36] 
S(t) 



v— In 

6 



4L 

— t sin 



( unt 
\2L 



(129) 



where v is the number of contact points and k v is a constant which depends on the 
geometry. This formula is in good agreement with the numerical data. For t <C L, it 
gives a logarithmic entropy growth. If L — > oo this persists for all times. 



Entanglement in solvable many-particle models 



41 



For numerical calculations, this is a more favourable situation. One can follow the 
evolution a much longer time. 

5.4- Periodic switching 

An interesting effect appears if one connects and disconnects two half-chains periodi- 
cally for a certain time r. One can call this a periodic local quench. Numerical results 
for the entanglement are shown in fig. 29. 




o ^ 1 1 1 1 1 1 1 

5 10 15 20 25 30 35 40 
t 



Figure 29. Entropy evolution for periodically connected chains and L = 40. Upper 
curve: r = 5, lower curve: r = 1. From [1]. Copyright IOP publishing, reprinted with 
permission. 

Features 

• Switching directly visible 

• Rapid switching: logarithmic increase 

• Slow switching: linear increase 

The curve for rapid switching resembles the result for a single quench, compare fig. 27. 
This can be understood as follows. The time-evolution operator for one period is 

U = U U 1 = e - iHoT e- iHlT (130) 

where H and Hi are the Hamiltonians for the the two configurations and do not 
commute. However, for small r, one can take the same Hamiltonian limit as for the 
transfer matrices in section 3 and combine the exponentials. Then 

U = e~ iR2r , if = i(if + #i) (131) 

The average time evolution therefore corresponds to a single local quench where the 
final system has a defect with reduced hopping amplitude t' = t/2 at the contact. For 
such a case, the evolution of S is similar as for a quench to a homogeneous system and 
the behavour is logarithmic in time. 

The curve for slow switching rises on average linearly. The interpretation is 
that here the disconnected system has enough time to "recover" and thereby the 
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entanglement gain repeats itself after each new connection. The problem can be treated 
analytically in a continuum model [37]. 

In general, one can express S in terms of the (time-dependent) cumulants of the 
probability distribution P n to transfer n particles. This provides a link to the so-called 
"full counting statistics" of the junction and thus in principle to measurable quantities. 
For the example given here, the distribution is Gaussian and only the second cumulant 
enters. 

5.5. Entanglement Hamiltonian and subsystem Hamiltonian 

The thermal form of the RDM automatically leads to the question whether there is a 
relation between % a and H a . In section 2.6 we have already seen that in general this is 
not so. But are there cases, where a relation exists ? 

The answer is yes. For example, it has been seen in Heisenberg ladders, where the 
subsystem was chosen as one of the two legs. We discuss here an example, which is 
somewhat simpler and a free-fermion model [38]. 
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Figure 30. Ladder geometry for a fermionic hopping model. 
The subsystem is chosen as one of the legs. 

Consider a hopping model on a ladder with opposite dispersion in both legs and 
hopping with amplitude 5 between them. The Hamiltonian is 

H = H 1 + H 2 + H' = J2l q a\a q - £ % b\b q + £ 5 (a% + b\a q ) (132) 



Diagonalizing (132) with a canonical transformation 



a q = u q a q + v q (3 q , b q = -v q a q + u q f3 q , u q + v q = 1, 



one obtains 



(133) 
(134) 



From that, one obtains the correlation matrix. Due to the translation invariance, it is 
diagonal in momentum space and in the subsystem 1 of the a's one has 



1 7 
( q =< a\a q >= v 2 q = -(l- 



This gives the single-particle eigenvalues 

'u q + % 



e q = In 



(135) 



(136) 
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and Tii has the form 



(137) 



If now the rung hopping 8 is large, one obtains e q = 2j q /5 and the relation 



o 



(138) 



This is a direct proportionality between the two Hamiltionians. 

Remarks 

• Holds for dominating rung couplings 

• Follows from first-order perturbation theory in Hi + H 2 

• Entanglement near maximum S = L In 2 

• Entropy extensive due to long interface 

For arbitrary S, the single-particle energies e q and 7 g are not proportional to each 
other. Therefore the hopping range in "Hi is in general different from that in H\. 

5.6. Concluding remarks 

I have given an account of the entanglement properties of solvable models, either free 
particle or integrable, and shown in particular that 

• One is lead to a thermodynamic problem 

• A particular Hamiltonian enters 

• Its spectrum determines the Schmidt weights 

• The ground states of homogeneous chains are weakly entangled 

• Global quenches lead to strongly entangled states 

Almost all considerations had to do with lattice models. These are the systems 
one studies in numerical investigations motivated by solid state physics or uses in 
quantum information. They also have the advantage that no divergencies appear in 
finite geometries. 

This does not mean that continuum systems are unimportant. The first calculations 
of entanglement entropies took place in the context of black-hole theory and thus in a 
continuum setting. And the use of conformal invariance has not only shown a deeper 
connection between the various models but also allowed to derive many special results. 
But that would be a lecture series in its own. Those who are interested can find a lot 
of material in a special issue of J. Phys. A 42 (2009). There, entanglement for free 
quantum fields is reviewed by Casini and Huerta [39] and within conformal field theory 
by Calabrese and Cardy [40]. 
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